home *** CD-ROM | disk | FTP | other *** search
/ Collection of Tools & Utilities / Collection of Tools and Utilities.iso / sound / fftscop4.zip / HILFE.TXT < prev    next >
Text File  |  1994-05-30  |  10KB  |  242 lines
  1. SMALLFFT.EXE   eine Abgemagerte Version um die Programmierung  
  2. SMALLFFT.C     besser verstehen zu können!
  3.  
  4. FFT.ZIP        die FFT-Routinen, alt (1988) aber gut!
  5.            jedoch mangelhaft dokumentiert...
  6.  
  7.  
  8.  
  9.  
  10.            Kurzbeschreibung FFT.ZIP
  11. ---------------------------------------------------------------------------
  12. die aufgenommenen Werte werden im real[] Array abgelegt,
  13. das imag[] Array wird mit Nullen aufgefüllt...
  14.  
  15. dann werden die Globalen Variablen samples und power gesetzt:
  16. samples muß ein Potenz von 2 sein: 2,4,8,16,32,64,128,256,512,1024,2048,...
  17.  
  18. int     samples, power;
  19. double  real[2048], imag[2048], max;
  20. samples = Anzahl aufgenommener Werte im real[]/imag[] Array;
  21. power =  log10((double)samples) / log10((double)2.0);
  22.  
  23. Jetzt wird die Funktion fft() ausgeführt!
  24. dieses ist sehr schnell, extrem viel schneller als die normale FT.
  25. das Real/Imaginär-Spektrum befindet sich danach im real[]/imag[] Array.
  26. Allerdings recht merkwürdig durcheinandergewürfelt (das ist bei FFT normal)!
  27. die permute(n) liefert uns die Position des Spektrums n in real[]/imag[]
  28. Beispiel: 
  29.               t=permute(n) 
  30.               r=real[t]
  31.               i=imag[t]
  32.               Betrag=sqrt(r*r + i*i) 
  33.               return Betrag;
  34.  
  35. dieses Beispiel ist die Funktion magnitude(n), welche den Betrag liefert!
  36. Wenn wir uns eh nur für das Betragspektrum interessieren benutzen wir
  37. also magnitude() und haben mit permute() nix mehr am Hut!
  38.  
  39.  
  40.  
  41. Iss auch ne Rückwärts (inverse) Fast-Fourier-Transfo möglich?
  42. Aber klaro, laut dem schlauen FFT-Transformation Buch das ich mal kurz
  43. in den Fingern hatte ist dazu lediglich der Imaginärteil zu invertieren, 
  44. (also mit -1 zu multiplizieren), fft() aufrufen, und dann wieder den
  45. Imaginärteil invertieren.
  46. aber dabei ist die Richtige Reihenfolge im Array wichtig (siehe permute)!
  47.  
  48. fft();   
  49. for(i=0;i<samples;i++) imag[i]*=(-1);
  50. fft();
  51. for(i=0;i<samples;i++) imag[i]*=(-1); 
  52.  
  53. nach diesem kurzen Programm stehen nicht etwa wieder die ursprünglichen
  54. Werte im Array, sondern Datenmüll, weil die Verwürfelung ignoriert wurde!
  55. Richtiger wäre:
  56.  
  57. fft();   
  58. for(i=0;i<samples;i++)
  59.               {
  60.                t=permute(i)
  61.                real2[i]=real[t];
  62.                imag2[i]=imag[t];
  63.               }
  64. for(i=0;i<samples;i++) real[i]=real2[i], imag[i]=imag2[i];
  65. for(i=0;i<samples) imag[i]*=(-1);
  66. fft();
  67. for(i=0;i<samples) imag[i]*=(-1); 
  68. for(i=0;i<samples;i++)
  69.               {
  70.                t=permute(i)
  71.                real2[i]=real[t];
  72.                imag2[i]=imag[t];
  73.               }
  74.  
  75. Nach der Ausführung dieses Programs stehen die urspünglichen Werte mit
  76. kleinen/grossen Fehlern wieder im real[]/imag[] Array und im 
  77. Zwischenspeicher-Array real2[]/imag2[] finden wir das Real/Imaginärspektrum!
  78.  
  79.  
  80. welches Format hat denn nu dasch Real/Imaginär/Betragsspektrum?
  81. 0 = Gleichanteil ....... samples-1 = samplerate
  82. leider können wir das Array nur bis zur Hälfte (samples/2 - 1) nutzen weil
  83. ein gewisser Nyquist das sogenannte Nyquist-Kriterium aufgestellt hat,
  84. nach dem bei Abtastung maximal Frequenzen der halben Abtastrate abgetastet
  85. werden können. 
  86. Erklärung mittels Systemtheorie: (bahh, würg, igittigitt)
  87. Ich versuche es unmathematisch visual anschaulich zu erklären,
  88. also Zeichblatt und Stift rauskramen und versuchen das Folgende zu zeichnen:
  89.  
  90. wir betrachten die Abtastung im Freqenzbereich:
  91. die kontinuierliche Abtastung (Kammfunktion im Frequenz- und Zeitbereich)
  92. mit Abstand   samplerate (Hz) im Freqenzbereich und
  93. mit Abstand 1/samplerate (s)  im Zeitbereich
  94. wird im 
  95. Freqenzbereich gefaltet mit dem Eingangssignalspektrum und im
  96. Zeitbereich multipliziert mit der Eingangsignal.
  97.  
  98. Faltung ist normalerweise besch... zu zeichnen aber wegen den Diracstößen
  99. und der Kammfunktion wirds viel einfacher, wir brauchen bloss das
  100. Eingangssignalspektrum um jeden Diracstoß zu zeichnen, als wenn jeder
  101. Diracstoß f=0Hz hätte, wenn wir das gezeichnet haben, 
  102. verstehen wir nicht nur das Nyquist-Kriterium, 
  103. wir wissen auch wie unser Spektrum zwischen samples/2 ... samples aussieht!
  104. Das sind nähmlich nicht weiter als die negativen Frequenzen des um
  105. Samplerate gefalteten Eingangssignalspektrums.
  106.  
  107. Also entspricht der linke Rand (0) und der rechte Rand (samples)  des Arrays 
  108. (permute beachten) dem Gleichanteil (0 Hz) des Eingangsspektrums, da der
  109. maximale adressierbare Wert im Array jedoch samples-1 ist finden wir den 
  110. Gleichanteil nur bei real[permute(0)],imag[permute(0)]... 
  111.  
  112.  
  113. Realteil:
  114. das Signal von 0...samples/2 - 1 wird um samples/2 gespiegelt.
  115. Imaginärteil:
  116. das Signal von 0...samples/2 - 1 wird um samples/2 gespiegelt, und zusätzlich
  117. invertiert (also das Vorzeichen gewechselt).
  118.  
  119. Betragspektrum: um das Betragsspektrum in den Zeitbereich zu transformieren
  120. zerlegen wir es in eine gerade und eine ungerade Funktion,
  121. der gerade Anteil kommt nach real[]
  122. der ungerade Anteil kommt nach imag[]
  123. dann wird Realteil und Imagiärteil wie angegeben gespiegelt!
  124.  
  125. für die inverse FFT ist natürlich noch die weiter oben erwähnte doppelte
  126. Invertierung des Imaginärteils notwendig!!!
  127. -----------------------------------------------------------------------------
  128.  
  129.  
  130. Informationstragend ist dabei eigendlich nur der Bereich von 0...samples/2-1
  131. entsprechend f=0 bis f=samplerate/2 also von 0 Hz bis zur halben Abtastfrequenz.
  132.  
  133. Wie meine Experimente mittels Soundkarte und Funktionsgenrator zeigten,
  134. ist beim Realspektrum und beim Imaginärspektrum das Vorzeichen dauernd am
  135. kippen, schade...
  136. das Betragsspektrum wird mittels Arctan(imag[permute(i)]/real[permute(i)]) 
  137. berechnet, leider eine ungerade und damit vorzeichensensitive Funktion...
  138.  
  139. das Betragsspektrum ist nicht nur sehr aussagekräftig, es ist auch wegen
  140. der Quadratbildung nicht vorzeichensensitiv also stabil!!!
  141.  
  142. es gibt da allerdings noch ein paar Probleme auf die ich gestoßen bin:
  143. bedingt durch die Abtastung kommt es wenn sich ein Vielfaches der 
  144. Abgetasteten Frequenz sich der Abtastrate nähert zu einer Schwebung
  145. (Amplitudenmodulation), daraus resultieren mit Sicherheit Fehler!
  146.  
  147. die FFT hat auch noch ein ähnlich gelagertes Problem:
  148. wie wir aus unserer Zeichnung ersehen können, ist die diskrete 
  149. Fourier Transformation bedingt durch die Abtastung 
  150. im Freqenzbereich periodisch, daraus ergibt sich der Umkehrschluß das
  151. auch irgendeine Art von Periodischen Verhalten im Zeitbereich vorliegt!
  152. leider trifft dieser Umkehrschluß zu !!!!!!!
  153. wenn irgendein Vielfaches der Periodenlänge der Abgetasteten Funktion sich
  154. der Zeit (samples/samplerate) also (Anzahl_Werte/Abtastfrequenz) nähert
  155. dann werden die Spektralinien hoch und schmal, wenn dies nicht zutrifft
  156. werden die Spektallinien entsprechend der Nichtübereinstimmung flach und
  157. Breit, die Energie liegt also in der Fläche!
  158.  
  159. Außerdem steht zu befürchten das sich die Beiden Probleme überlagern
  160. (belagern) weil sie nähmlich unter den selben Bedingungen auftreten.
  161. Möglicherweise handelt es sich hierbei nicht um eine Überlagerung zweier
  162. Probleme sondern lediglich um ein von Zeitbereich in den Frequenzbereich 
  163. transformiertes Problem, verwandt sind bei Probleme auf jeden Fall!
  164.  
  165.  
  166.  
  167.  
  168. Es gibt zwei möglich Lösungen dieses Ärgerlichen Problems,
  169.  
  170. 1. man benutzt eine geeignete Fensterfunktion windows(), daß heißt:  
  171.    man multipliziert das Signal im Zeitbereich mit einer sanft Anklingenden
  172.    und sanft ausklingenden Funktion (sin,gauss...) so das zum Beispiel 
  173.    eine Abgetastete Funktion nicht schlagartig anfängt und abbricht
  174.    sondern langsam ein und wieder ausgeblendet wird.
  175.  
  176. 2. man nimmt das abzutastende Signal mehrfach mit unterschiedlichen 
  177.    Abtastraten auf und überlagert die entstehenden Betragsspektren!
  178.    dabei empfehle ich instinktiv das Quadratische Mittel zur Überlagerung
  179.    heranzuziehen. Natürlich müssen die Frequenzen umgerechnet werden,
  180.    und genau da liegt ein weiterer gefährlicher Fallstrick, denn
  181.    die Mega-Blaster 16 (Mozart) Soundkarte emuliert ja nur eine 
  182.    Soundblaster Pro Soundkarte, das heisst das sich aufgrund der 
  183.    unterschiedlichen Hardware und des Begrenzten Wertebereich (8.Bit) 
  184.    der Timerconstante kleine Unterschiede betreffend die tatsächliche 
  185.    Abtastfreqenz ergeben. Die Folgen für ein Überlagerung wären Katast